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We analyze strongly interacting Fermi gases in the unitary regime by considering the generalization to an 
arbitrary number A'^ of spin- 1/2 fermion flavors with Sp(2N) symmetry. For TV ^ oo this problem is exactly 
solved by the BCS-BEC mean-field theory, with corrections small in the parameter 1 /N. The large- expansion 
provides a systematic way to determine corrections to mean-field predictions, allowing the calculation of a 
variety of thermodynamic quantities at (and in proximity to) unitarity, including the energy, the pairing gap, and 
the upper-critical polarization (in the case of a polarized gas) for the normal to superfluid instability. For the 
physical case of TV = 1, among other quantities, we predict in the unitarity regime, the energy of the gas to be 
^ — 0.28 times that for the non-interacting gas and the pairing gap to be 0.52 times the Fermi energy. 
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I. INTRODUCTION 
A. Motivation and Background 

Recent advances in atomic gases near a Feshbach resonance 
(FR) have led to the experimental realization of resonantly- 
paired superfluidity of fermionic atomic gases''^'"'''*'^'^. Al- 
though atomic interactions are generally weak, such superflu- 
idity has been achieved using the strong attractive interactions 
provided by proximity to a magnetic-field tuned Feshbach res- 
onance between the two hyperfine levels (isomorphic to spin- 
1/2 states) undergoing pairing. Indeed, the tunability of such 
Feshbach resonances, experimentally controllable via an ex- 
ternal magnetic field, allows unprecedented access to a wide 
range of fermionic interaction strengths characterized by the 
vacuum s-wave scattering length Os that diverges for a reso- 
nance tuned to zero energy. As the FR detuning S cx — 1 /a^ is 
varied, the character of superfluidity evolves from a weakly- 
paired Bardeen-Cooper-Schrieffer(BCS) regime at large pos- 
itive FR detuning (where < 0) to a strongly-paired molec- 
ular Bose-Einstein Condensate (BEC) regime at large nega- 
tive detuning (where as > 0). The experimental signatures 
of such pairing and superfluidity in ultracold gases of ''°K and 
^Li have included direct measurements of the condensate den- 
sity'-^ and the pairing gap^ and the observation of vortex lat- 
tices in rotating clouds^. 

From a theoretical perspective, a quantitative descrip- 
tion of such resonantly-interacting Fermi superfluids is well- 
developed away from the resonance, where is short com- 
pared to atom spacing. Then, a controlled perturbative expan- 
sion in a natural small parameter |as|n^/"^ ^ 1 (with n the 
fermion density) allows a quantitative theoretical analysis of 
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regimes. 



How- 



both the deep BCS^ "' and deep BEC 
ever, the aforementioned present-day experiments are typi- 
cally in the crossover between the BCS and BEC regimes, 
where |a.,|n^/'^ ^ 1. Thus, the absence of a small parame- 



j.gpi4ji5ji6 jjg^j. jjjg resonance precludes a systematic perturba- 
tive expansion for a quantitative description of this theoreti- 
cally challenging regime. 

A particular point of interest is the so-called unitarity point, 
precisely at zero FR detuning, at which aj^ =0 and fermion 
scattering is characterized by a maximum scattering phase 



shift of it/ 2. At this point, the system does not contain any 
scale besides the Fermi wavelength kp^ set by the atom den- 
sity n — kp/{37T^), and the only energy scale is the Fermi 
energy ep = h'^kp/2m and therefore the free energy, and any 
quantities related to it, is given hy ep multiplied by a univer- 
sal dimensionless function of kBT/cp—-. In particular, this 
implies that, at zero temperature, the internal energy per parti- 
cle e is simply proportional to that of a non-interacting Fermi 
gas 



(1.1) 



Recently, there has been much theoretical interest in comput- 
ing ^ and other universal parameters, motivated by the pos- 
sibility of attaining a quantitative understanding of unitary 
quantum gases as realized in cold-atom experiments but also 
having applications to nuclear physics and astrophysical sys- 
tems such as neutron stai-s' ^i''^i^°i^'i^^i^^'^'^'^^i^^-^ ^-^^-^'^. The pur- 
pose of this paper is to provide a theoretical framework for a 
systematic determination of the value of such universal param- 
eters at unitarity. Our method is based on the introduction of 
an artificial small parameter 1/TV, with TV the number of dis- 
tinct "spin"- 1/2 fermion flavors with Sp{2N) symmetry. We 
then use it to extrapolate to the experimentally-relevant case 
of a single flavor (TV = 1) of two opposite-"spin"-l/2 (hyper- 
fine levels) fermionic atoms. 

The motivation for such a generalization is that, for TV — > 
00, the problem may be solved exactly-^", with the solution 
taking the form of the standard BCS mean-field theory. We 
can then compute corrections in the small parameter 1 /TV, ob- 
taining a systematic expansion (in principle to arbitrary order) 
about this solvable limit. Since it is believed that there are 
no phase transitions with decreasing TV (the large-TV solution 
having the same broken symmetry as the exact ground-state 
at TV = 1), we expect the large-TV results to be smoothly con- 
nected to the physical case of TV = 1. The large-TV expansion 
may be thought of as a way to systematically organize correc- 
tions around the well-known BCS mean-field solution. Since 
this TV — > CX3 solution is known to give a relatively good esti- 
mate, we expect the 1 /TV-expansion to converge rapidly. 

To further motivate our study, it is instructive to briefly re- 
view other theoretical approaches. As we have noted, near 
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the unitarity point a quantitative theoretical understanding of 
this strongly correlated problem is limited by the absence 
of any physical small parameter-i^>i^. However, since, the 
thermodynamic quantities are expected to evolve smoothly 
from the BCS to the BEC limits, a qualitative description 
has come from various uncontrolled (by any small parameter) 
schemes, that interpolate between these two regimes. Start- 
ing from a variational wavefunction (equivalent to mean-field 
theory). Eagles'^ and Leggett'" studied the crossover from 
the BCS to the BEC limit at zero temperature. Nozieres 
and Schmitt-Rink studied the same problem at finite tem- 
perature, taking into account the attraction of fermions be- 
yond mean-field theory^' . The repulsion between bound pairs 
was first worked out using a functional integral formalism.^^-^^ 
and a self-consistent theory^. However, except for numeri- 
cal Monte Carlo calculationsi2i22i2iiS?., a full description of the 
system around unitarity in a quantitative fashion has remained 
elusive due to the absence of a small parameter. 

More recently, stimulated by progress made in the theory 
of critical phenomena (where, similarly, no physical small pa- 
rameter exists near a critical point)"*^, progress has been made 
by studying the system in d dimensions, treating the unitary 
regime in a systematic expansion in a small parameter associ- 
ated with the dimensionality of spaca2^i212^i2^. 

The large-A^ approach used here is close in spirit, but is 
complementary to such an expansion-in-dimension study, in 
that we also introduce an artificial, but distinct small param- 
eter and perform a systematic perturbative expansion in it. 
Such an approach has also been extremely successfully in a 
variety of field theory and statistical physics contexts^il^i^, 
applied to a description of continuous phase transitions close 
to a critical point. 



B. Summary of results 

Our results are predictions for the lowest nontrivial order 
in for an array of physical quantities at and around the 
unitary point. At unitarity aj^ ~ 0, a Fermi gas at den- 
sity n per fermionic flavor is characterized by a single energy 
scale ep — {STT^n)"^^^ /{2m). In this regime, we find that, in 
the symmetry-broken (paired-superfluid) phase at T = 0, the 
chemical potential /i, the order parameter A, and the excita- 
tion gap Aoxc to order 1 /N are given by 



A/ep = 0.6864 - 0.163/iV + C'(l/iV2), (1.2) 
AcxcAf = 0.6864 - 0.196/iV + 0{l/N^), (1.3) 
Ai/eF = 0.5906 - 0.312/iV + ©(l/iV^), (1.4) 



tuting = 1 in Eqs. il.2\ and ( ll.4l l gives our predictions 



Aoxc/ 



0.523, 
0.490, 
0.279. 



(1.5) 
(1.6) 
(1.7) 



We have also computed the critical temperature Tc and the 
chemical potential /i at Tc to be given by 



ksTc/ep 
li/ep 



0.4964 
0.7469 



1.31/iV + 0(l/iV2), (1.8) 
0.58/iV + O(l/iV2). (1.9) 



The 1 /A^-expansion can also be applied to the problem of 
a unitary Fermi gas with population imbalance 6n — ~ 
between the densities rij, rii of the two hyperfine fermion 
components, that has been realized in experiments^SilLi^ and 
has been a focus of intense theoretical activity; see e.g., Refl43l 
and references therein. Mean-field theory predicts^^iimi^, 
that, for arbitrarily small polarization P = 5n/n, a Fermi 
gas in the unitary and positive-detuning regimes phase sepa- 
rates into a polarized normal (unpaired) phase and an unpolar- 
ized paired-superfluid phase, with the average density and po- 
larization equal to the experimentally imposed values. Such 
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with the first N ^ oo term corresponding to the well-known 
mean-field theory result?.. As we discuss in the Appendix, 
although A = Acxc in the N — oo limit, at subleading order 
these quantities differ, as seen in Eq. ( 11.21 ). 

By using scaling arguments, it can be shown that, at unitar- 
ity, the ratio iJ,/ep is equal to ^, defined by Eq. ( II. lb . Substi- 



FIG. 1 : Chemical potential fi and order parameter A as a function of 
(kpas)^^. The dashed lines are the N ^ cx) results for /i and A. 



The diamond symbols include the 0{1/N) corrections, evaluated at 
N — 1. The solid lines are a guide to the eye. The star symbols 
at unitarity are the results of quantum Monte Carlo calculation from 
Refill 
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phase separation has, in fact, already been observed experi- 
mentally. We find that, within the 1/A^-expansion (consistent 
with mean-field theory) at unitarity, the paired superfluid is 
unstable to phase separation for an arbitrarily small polariza- 
tion, i.e.. Pel = 0."^^ 

For sufficiently large P, or, equivalently, sufficiently large 
applied "Zeeman field" h — -^(/^i — fii) a uniform polarized 
normal phase is stable and the corresponding critical values 
Pc2 and hc2 are also universal quantities at unitarity. We find 
their 1/iV-expansion to be given by 

hc2/eF = 0.6929 + 0.087/iV + 0{1/N^), (1.10) 
Pc2 ^ 0M26 - Omi/N + 0{1/N^), (1.11) 
fi/ep = 0.8585 - 0.458/iV + 0{1/N^), (1.12) 

where the chemical potential is evaluated at hc2- 

This paper is organized as follows. In Sec. |II] we introduce 
the one-channel model and its generalization to the N flavor 
case. In Sec. |III1 we show that the solution for cx) is 

equivalent to the mean-field solution of the iV = 1 model. 
We then take into account, systematically, the fluctuation cor- 
rections to leading order in 1 /N to the internal energy, exci- 
tation gap, the critical temperature and the upper-critical po- 
larization, and extend our results to the vicinity of the unitary 
point. In Sec. IIVI we discuss our results and compare them 
to recent measurements and predictions of other theoretical 
approaches. In Sec. [Vl we provide a brief summary. 



II. MODEL 

We study a three-dimensional resonant Fermi gas confined 
to a box of size L x L y. L ~ V , described by a Hamiltonianii 
(setting h^l) 



A / d3rV>j(r)V'|(r)V'i(r)^TW 



(2.1) 



where A < is the fermionic interaction and -0^ (r) , ip^ (r) 
are, respectively, the fermion creation and annihilation oper- 
ators at position r and hyperfine state ("spin "-1/2) a, which 

obey the usual anticommutation relation |?/'CT(r), V'!' — 

5{y ~ Y')5a.<j' ■ The chemical potential fjia fixes the average 
density of particles with spin a. For simplicity we now set 
/i| = III = fi, although later we shall allow the possibility 
of a nonzero chemical potential difference in order to study a 
polarized Fermi gas. 

The partition function Z = Tr [e^''^] gives the free- 
energy density 



/ = In Z, 



(2.2) 



where (3=1 /T is the inverse temperature, taking fcs = 1 
throughout. 



The two-particle scattering length is related to the 
strength of the atomic interaction A via the relation 



TO 



1 1 V ^ 

k 



(2.3) 



where ek = fc^/(2m) is the free fermion dispersion (kinetic 
energy), taking h = 1 throughout. In the second term on 
the right side of Eq. (12.3b an ultraviolet cutoff A ^ 1/ro is 
implied, set by the effective interaction range ro <C 1 /fcir'"* '^. 

We now proceed with the generalization of Eq. ( 12.11 ) to a 
model of N fermion flavors: 

N 



N 



+ /d^r^,\(r)0t^(r)V',i(r)0,TW, (2.4) 

that is clearly equal to Eq. (12. Il l for the case of iV = 1 and pos- 
sessing invariance under the symplectic group Sp{2N) (see 
Ref. 47). As discussed in Sec. II Al the benefit of this expan- 
sion to N fermion flavors is that an exact solution may be 
obtained in the large- limit. To see this, we write the par- 
tition function in terms of an imaginary-time coherent-state 
functional integral 



Dipi^ {x)Dilji„ {x) exp(-5). 



(2.5) 



with t|J^^{x),'^|Jicr{x) Grassmann fields labeling the corre- 
sponding coherent state at a space-time point x = (r, t), and 
with the imaginary-time action S given by 

+ dT I (f'T',p^^{x)^^i{x)ll)jl{x)^pjl{x). 



(2.6) 



We decouple the atomic interaction via a standard Hubbard- 
Stratonovich transformation, that relies on the Gaussian inte- 
gral 



e " 



c/>'(x)<p{x) 



— db{x)db*(x)exp(—b*{x)b{x) 
ttX J A 

+ b{x)<j)*ix)+b*ix)(t>ix)), (2.7) 



which we shall utilize with 0(.t) given by a superposition of 
N flavors of bilinear fermionic fields 



N 



(2.8) 



With this transformation, the fermion fields appear quadrati- 
cally and can therefore be formally integrated out. The par- 
tition function is then given by a functional integral over the 
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boson field b{x), 

Z = Z^^ j Db{x)Db* (a;)e-^[''l , (2.9) 



with the effective bosonic action given by 

S[b] 



A 

-NTi-ln{-G-^ [b{x),b*{x)]), (2.10) 

where the trace is over space, imaginary time, hyperfine and 
flavor states, the Green's function (written in Nambu repre- 
sentation for the fermions) is 



G-\x) 




b{x) 

2m 



(2.11) 



and the Hubbard-Stratonovich normalization factor Zb = 

J Db{x)Db* (x) exp (iV dr J d^rb* {x)b{x) /a) . Note 

that in the limit where TV = 1 the effective action reduces 
to the one-channel model derived by Sa de Melo, et al?^. 

We use this formulation to compute a number of thermody- 
namic quantities that characterize the resonant unitary Fermi 
gas. To this end, we observe that in the limit N oo, these 
expectations values are dominated by the saddle point of the 
functional integral. This is the origin of the aforementioned 
claim of an exact solution in the limit of large- A^. As is well- 
known, such a saddle-point evaluation amounts to the stan- 
dard BEC-BCS mean-field approximation to the one-channel 
model"^ that is only approximately valid for the case iV = 1. 
At large N , fluctuations around this saddle point are small in 
the parameter 1/iV, generating 0{1/N) corrections, that or- 
ganize into a systematic expansion in powers of 1/iV for the 
partition function and physical quantities derived from it. 

By representing the bosonic Green's function by a line, the 
partition function is given by a series of closed loops. The 
higher order terms in 1 /N can be classified according to the 
number of loops for a particular diagram. Thus the 1/N- 
expansion is equivalent to a so-called loop expansion. 

The implementation of the loop expansion is performed 
around the (possibly complex) saddle point A, where A is 
the mean-field value of b{x). At low temperature, that is our 
focus here, the saddle point yields a non-zero value for A. 
The complex field fluctuations b{x) around the saddle point 
are defined as 



-^b(x) = b(x) - A. 



(2.12) 



The Green's function matrix can be formally separated into 
the saddle point and fluctuation contributions, 

1 



G-\x) = G^J^ix) + -^G^^^ix), (2.13) 



where 



-dr 



A* 



' 2m 



(2.14) 



and 



GvMx) = - " ^ ' , (2.15) 
The bosonic action S[b] can be expanded in powers of 1 /N 

-Tr 

1 Ab*{x) + A*b{x) 



N ^ X 



dr / d-'r 





m— 1 



1 



m 



]\[m/2 



N 



Tr 



A 



G^o)Gl,][b{x)] 



N A 



(2.16) 



Up to this point, all the transformations of the partition 
function are exact. An approximation is made when only a 
finite number of terms in the action are considered, and the 
small parameter 1 /N provides a systematic way to organize 
this expansion. To analyze the problem to lowest nontrivial 
order in it is sufficient to include terms up to m = 2 in 
the action, yielding an approximate bosonic action 



N 



]_oil/N) 



where 
and 

g(l/N) ^ 



lAI 



Tr 



In 



(2.17) 



(2.18) 



dr / d-^r 





iTr 
2 



\b{x)\' 



Gl,][b{x)]G<^o)Gl,][b{x)]G^o) 



(2.19) 



that we shall analyze in the remainder of the manuscript. We 
note that, as usual, terms linear in b (the to = 1 terms) au- 
tomatically vanish by a virtue of the definition of the saddle- 
point as the extremum of the action S^^\ 



in. LARGE-TV EXPANSION 

In the present section, we use the results of Sec. |ll]to con- 
struct an expansion of the free-energy density / to leading 
order in small 1/TV of the form 



(3.1) 



where the leading-order term f^^'' comes from a saddle-point 
approximation to the bosonic functional-integral expression, 
Eq. (|2.9l l, for the partition function. 



A. Saddle-point approximation: Mean-field level 

In this section, we compute the lowest order approximation 
to the free-energy density, namely the first term of Eq. (13. 11 1. 
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that is exact in the N oo hmit. This corresponds to an 
evaluation of the free-energy density within the saddle-point 
approximation /''^^ = — -gW InZ'"' = ■JyS'--^\ which gives 



N ^ oo limit is given by 



^(0)^_, 



lAI- 



(£^k-a)- 



2 



(3.2) 

where we included an overall constant in the second term 
so that for A ^ 0, reduces to the free-en ergy density 
of a free Fermi gas. Here, £'k — + | Ap is the spec- 
trum of the quasiparticles and = Ck ~ A*. Notice that this 
N ^ oo limit yields a result for the free energy (and therefore 
all subsequent thermodynamics quantities) that is identical to 
the usual mean-field treatment based on a BCS ansatz for the 
ground-state wavefunction. 

The gap equation, determined by minimizing the free en- 



ergy with respect to the gap. 



a/ 

OA 



0, yields 



mA 



A ^ ^tanh(/3£;k/2) 1 



2ek 



7^ • (3.3) 



The total atom density is n x A^, where n is the atom density 
of each of the fermion flavors. It is determined through a La- 
grange multiplier via nN — which yields for — > oo 

n = i (^1 - tanh(/3£;k/2)) . (3.4) 

This number equation, along with Eq. ( 13.3b . describes a reso- 
nant Fermi gas within the mean-field approximation. 

The solution to these equations at zero and finite temper- 
ature has been discussed at length in, e.g., Refs. [l(illl4l3ll 
l33ll48[ For completeness, below we review this mean-field 
solution near unitarity at zero temperature, near the transi- 
tion temperature Tc, and its generalization to an imposed 
polarization^. 



1. Zero-temperature limit 

At zero temperature, the mean-field ground state at unitar- 
ity, [kpas)^^ = 0, is a superfluid characterized by the BCS 
wavefunction, with the gap and chemical potential given by 



/sS^^/ep = 0.6864, 
/i(°VeF = 0.5906, 



(3.5) 
(3.6) 



simply obtained by solving Eqs. (13.3b and ( 13.4b at T = 0. 
Here, the superscript and subscript refer to the 1 /N = and 
T — limits, respectively. 



2. NearT^ 

At finite temperatures, a second-order superfluid-to-normal 
transition takes place at a critical temperature Tc, that in the 



fcsri"VeF = 0.4965, 
^(o)/e^ = 0.7469, 



(3.7) 
(3.8) 



ji^c^ is the 



obtained by solving these equations at A = 0^ 
chemical potential at the transition as determined by the fixed 
density condition Eq. (13.41) . 



3. Polarized Fermi gas 

We can also study a polarized Fermi gas (one with popula- 
tion imbalance between the spin components) in the saddle- 
point approximation, achieved through coupling the fermions 
to an applied "Zeeman" field. Physically, this corresponds to 
imposing different chemical potentials, jiyi — iizt h for the 
up and down spins, and yields the free-energy density 



k 



tT=± k 



(3.9) 



The detailed mean-field phase diagram found by minimizing 
/(°) has been discussed in detail elsewhere'*^ '*^. The relevant 
parameters characterizing the phase diagrams are the detun- 
ing and h, or equivalently (and more directly related to exper- 
iments at fixed species imbalance) as a function of an imposed 
polarization 



df 
dh 

9/ 



(3.10) 



with the density imbalance m and density n. We are partic- 
ularly interested in the region near unitarity, where for suf- 
ficiently large h at fixed /x there is a first-order phase tran- 
sition between the superfluid (s) and normal {n) states. At 
fixed imposed density n, the first-order transition opens up 
into a regime hd < h < hc2 of phase separation, where 
on the BCS side and around the unitary point a gapped un- 
polarized superfluid and a polarized Fermi gas coexist. For 
h < hci, the system is in the homogeneous fully-paired su- 
perfluid state. At unitarity, the lower-critical polarization for 
entering the regime of phase separation is P^i = For 
h > hc2, the normal state is stable. The instability to phase- 
separation from the normal state takes place at h — hc2 when 
the free energies of the normal in) and superfluid (s) states 
cross, subject to the constraint on /i that the normal-state num- 
ber equation is satisfied. Thus, we equate fs = fn, using the 
N ^ oo approximation to the superfluid-state free-energy 
density fs ~ /'"'(A) (with A determined by the minimiza- 
tion dfs/dA = 0), and the normal state free-energy density 
given by /„ — /s(0). Solving these equations at zero temper- 
ature, we find fi, h and A at the transition (in agreement with 
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earlier wor lr^' ° ) 



^l'^^/eF = 0.8586, 
^SV^F = 0.6930, 
aSV^f = 0.9979. 



(3.11) 



The upper-critical polarization is given by combining the 
normal-state polarization Eq. (13.10b with / given by /„ (tak- 
ing T = 0) 



P 



(3.12) 



{i + h/iif/^ + {i-h/iiY/^' 

with the location of the transition defined by Eq. (13.111 1. yield- 
ing the mean-field result 



Pc2 = 0.93261, 



(3.13) 



at unitarity. 

The quantities that we have computed in this section within 
the iV — > oo limit (that is equivalent to the BCS mean-field 
theory), are universal, i.e., independent of the microscopic in- 
teractions and can be obtained as derivatives of a universal free 
energy that is given by the Fermi energy (set by the atomic 
density n) times a system-independent scaling function of di- 
mensionless variables such as T/tp and h/ep- For example, 
taking h — for simplicity, based on dimensional grounds, 
the free-energy density at unitarity can only depend on T and 
the density n (which we represent through €p, see Eq. dl.ll ) ) 



kBT 



(37r2)n2/3/2m 



(3.14) 



where g{x) is a dimensionless function. At zero tempera- 
ture, the internal energy per particle [quoted in Eq. (Il.lb l 
e — f{T — 0)/n is ^g{0)ep, while the chemical potential 

/i = 1^ is equal to /i = g{0)eF- This leads to e = 
that together with Eq. ( 13.6b gives our lowest order — > oo 
(mean-field) estimate of the parameter ^ : 



^ = £,(0) « 0.5906. 



(3.15) 



Other quantities that are connected to ^ at T = are the bulk 
modulus B = n^d^f/dn^ given by: 



(3.16) 



and the first- (isothermal) sound v = ^B/nm. 

In a similar fashion, one can derive the equation of state 
of the system and determine relations between the entropy, 
pressure and internal energy"^^. Since these relations follow 
from scaling arguments, they hold for any number of fermion 
flavors. Consequently, these relations are expected to be pre- 
served order by order in the 1/A^-expansion. 

Having obtained leading-order results in 1 /N (already re- 
ported in the literature as they correspond to the standard 
mean-field theory, exact for the N ~ oo model), in the next 
subsection we proceed to calculate the leading-order correc- 
tions in small 1 /N. We shall see that the large- expansion 
of the one-channel model organizes subleading corrections to 
mean-field theory, essentially amounting to the random-phase 
approximation. 



B. Leading order in 1/A: Random Pliase Approximation 
(RPA) 

We now consider the leading-order corrections in 1/N to 
the results of the previous section, which requires an evalu- 
ation of the subleading contribution to the free energy. This 
is given by a Gaussian integration over b fluctuations around 
the saddle point, governed by the effective action Eq. ( 12. 19b . 
explicitly given by 

^ - 2 ^ ^^~'^^> \B*{q) A{-q)) \b*{-q)) ' 

(3.17) 

where q — (q, il^) and the bosonic Matsubara frequency 
ilf = 21x11 fi with I an integer The matrix elements of the 
polarization matrix are 



1 



1 



A ' ^^E(G(o)(fc + '?))ii(G(o)(fc))22' 



1 



(3.18) 
(3.19) 



which satisfy the relations A{q) — A*{—q), B{q) — B{~q). 

The quantity G(o) appearing in Eqs. ( 13.18b and ( 13.191 ) is 
the standard saddle-point approximation to the single-particle 
Green's function in the BCS state 

~1 ( iuje + -A 
-A* 



G(o)(fc) 




2 \ 
k V 



iuJjj—Ek 



iuje+Ek 



iujf — Ek 



iLj£-\-Ek iujf—Ek 



iuj,+Ek 



where 




uu = \ - I I 



(3.20) 



(3.21) 



(3.22) 



are the usual BCS coherence factors and k = (k, ioje) and 
uJe = {2£ + l)7r//3 are the fermionic Matsubara frequencies 
and for simplicity we have taken A to be real. Evaluation of 
these normal and anomalous particle-particle bubbles yields, 
respectively, 

Mq) = -T + ^y, (3-23) 



A 2V 



(tanh(/3£;k/2) + tanh(/3£;k+q/2)) 



^k+q"k 



9 9 
<+q^k 



\ ^-k+q + -Bk — iflg £^k+q + -Bk + i^e 
+ (tanh(/?£:k/2) - ta.iih{l3Ek+j2)) 



"k+q"k 



"k+q^k 



Eu 



^k+q 



7 



2V 

K 

(tanh(/3£;k/2) + tanh(/3£;k+q/2)) 
1 1 



(3.24) and n 



SA 



Ek+q + Ek + i^e Ek+q + Ek — i^t 
(tanh(/3£:k/2) - tanh(/3£;k+q/2)) 
-1 -1 



Ek + -Ek+q - i^^e -Ek + Ek+d + ii^t 



In terms of these the Gaussian functional integral over the 
b{x) and b* (x) fields finally gives the 1 /N contribution to the 
free-energy density 



/ 



(l/JV) _ 



1 



213V ^ 

9 

(3.25) 

where the function is evaluated at an imaginary time (5 = 0^. 

When inserted into Eqs. ( 13.1b . and using Eq. ( 13.21 ). 
Eq. ( |3.25t yields an explicit expression for the free-energy 
density of a resonantly-interacting Fermi gas to leading order 
in 1 /N . We now proceed to compute the effect of the 1 /N 
corrections in various limiting regimes, starting with T = 0. 



(0) 
Mo . 



Solving for 5A and 5^ yields 



(0) 
(0) 



5^^(i/iV) ) ' (3-31) 



thus expressing the shifts and in the chemical potential 
and gap in terms of various derivatives of the integrals e'"' 
and e^^/^-'. At unitarity, these can be calculated numerically 
and are given by 



a^^e(o) = -1.0804n/eF, 
S^AE^"^ -1.2556n/eF, 

Saae^"^ = 1.0804n/eF, 
a^£(i/JV) ^ _o.542n, 
9^£(i/JV) ^ _o.216n. 



(3.32) 
(3.33) 
(3.34) 
(3.35) 
(3.36) 



Substituting Eq. ( |3.31| l and the appropriate parameters into 
Eqs. ( 13.281 ) and ( 13.291 ) gives the 1/N corrections to /i and A 
quoted in Eqs. ( I1.2l i and ( 11.41 ) of the Introduction. 



2. Near Tc 



1. Zero-Temperature Limit 

At zero temperature, the energy density (to subleading or- 
der) is 



- = + -e 

N N 



0) and 



(l/JV) 



(3.26) 



where = f°\T 



1 



Next, we consider the vicinity of the transition temperature 
Tc, for which the saddle-point {N oo) results were pre- 
sented in Sec. lIII A"2] The 1/iV corrections to the critical tem- 
perature can be evaluated by computing the Thouless criterion 
(c'a/ — 0^) and particle number equation from the free en- 
ergy. Combining the two equations allows the determination 
of the critical temperature and chemical potential correction. 
Writing the saddle-point solution with the indices \ we can 
parameterize the 1 /N corrections to T and /i near Tc via 

T = ri°) + ^JT + . . . , 



(3.27) 

At large N, the chemical potential ji and gap parameter A 



(0) 



N 
1 

TV 



(5/i - 



(3.37) 
(3.38) 



are very close to their N 
order 0(1/^^) : 



oo values, within a correction of 



A(°) 



- — (5A 

N 



(3.28) 
(3.29) 



analogously to the zero-temperature case in the preceding sub- 
section. As in that case, the Thouless criterion and particle 
equation yield a set of two coupled equations 



5T 



a^/(i/JV) 



(3.39) 



Inserting Eqs. (13.281 1 and ( 13. 29b into Eq. ( 13.271 ) and using the 
gap = de/dA and number nN = —de/dfj, equations, we 
find that 6 A and Sfi satisfy 



where all the derivatives are evaluated at T = Ti"^ and /i 
fi^^ and A = 0+. Solving for 6T and Sfi, we obtain 



(0) 



9au£(°) aAAe(°W UA 



6fi 



(3.30) 



ST 



where all the derivatives (indicated by the shorthand notation 



dab 



d_d_ 

da db 



) are evaluated at the saddle point, i.e. A 



(3.40) 
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At unitarity, we find the matrix elements of Eq. ( I3.40l i to be 



f(0) 
J 


= -0.678 ncF, 


(3 41) 




= —1 1720n/ep 


G 42) 




— —1 2581r?/f F 


43) 




1 f;fiQ1 ri A /*:^ 

— — 1 .ouy 1 /tz_i/ 1^, 






= 2.3608?iA/e^, 


(3.45) 




= -O.TOlneF, 


(3.46) 




= -2.339n, 


(3.47) 


a^/d/iv) 


= 2.195nA/ei.. 


(3.48) 



Using Eq. ( I3.40l i together with these matrix elements inside 
Eqs. ( 13.37b and (13.38b gives the 1/iV correction to the critical 
temperature and the chemical potential quoted in Eqs. ( |1.8b 
and ( |1.9b of the Introduction. 



3. Polarized Fermi gas 

We now determine the 1 /N corrections, beyond the mean- 
field result reviewed in Sec. 1111 A 31 to the upper critical po- 
larization Pc2- Recall the two criteria to determine the upper- 
critical polarization, above which the normal state is stable. 
These are given by equality of the superfluid and normal- 
state free energies (the first-order transition condition) and the 
normal-state number equation, that are, respectively, given by 



fniP-.h) 



Nn 



fsili.KA), 



(3.49) 
(3.50) 

I 



Together, these determine the chemical potential difference 
below which the normal state is unstable to phase separation. 
In our formulation of the theory, the superfluid-state energy 
on the right side of Eq. ( 13.491 ) is a minimum with respect to a 
variational parameter, namely A, 



dA 



0. 



(3.51) 



Therefore, there are three parameters to solve for, namely 
A, /I, and h, using these three conditions. Once these have 
been determined, the upper critical polarization follows from 
m/n, where m is the density imbalance. In 



Pc. 

the superfluid state, the single particle excitations are gapped 
(for h < A). As a consequence, at zero temperature, the sus- 
ceptibility to the polarization field h vanishes, therefore, the 
energy of the ground state is independent of h and is given 
by Eq. ( 13.261 ). In the case of the normal state, the result is 
also relatively simple as the polarization field corresponds to 
a mere shift of the chemical potential. The correction to 
the action is 



^(i/^)=^S*(g)r-i(g)6(g), (3.52) 



where 



r-'(<z) = - 



A AV 



E 



Ea=± tanh(/3(gi,+q/2 + ah)/2) + tanh(/3(gk-q/2 + <jh)/2) 

- ^k+q/2 ^ ^k-q/2 



(3.53) 



The free-energy density contribution due to fluctuations is 
given by 



J n 



(l/JV) 



i-^e^"^^ln[-Ar-i(g)] . (3.54) 



The zeroth order /l"^ result can be obtained by taking the limit 
A ^ in Eq. ([33 

fn^ = E E (^k - I In [2 cosh(/3(a + <jh)/2)]^ . 

^ (3.55) 
As in the preceding subsections, the parameters /i h, and 
A acquire corrections, due to Gaussian fluctuations, that are 



small in and can therefore be written as 



A ^ 


^c2 






Pc2 ^ 




h = 


"'c2 ^ 





(3.56) 
(3.57) 
(3.58) 



Inserting this parameterization into Eqs. ( 13.491 ), ( 13.50b . and 
(13.51b . we find a set of three coupled equations giving these 
corrections 



9 



fd,i"^-dj'^ aA/r-9A/. 



(0) 
(0) 



dAhf, 



(0) 
(0) 



9f,Afi 

9a a/; 



(0) 
(0) 



Sh 



Al/N) 
' n 

dAf. 



Jl/N)\ 



(1/N) 

1 

(1/N) 



(3.59) 



where all the derivatives are evaluated at A = Aj,2 , h = ' 
and /i = fif^ ■ Restricting attention to zero temperature, many 
of the matrix elements vanish because : (1) at unitarity, the 
superfluid state is gapped and its susceptibility to a "Zeeman" 

I 



The matrix elements, at unitarity, are 

9^ei°) = -1.7527n, (3.64) 

d^^Aef'^ = -1.5140n/eF, (3.65) 

Oaab'-P = l.S027n/eF, (3.66) 

a^£l°) = -n, (3.67) 

dhe'^^ ^ -0.9326n, (3.68) 

5mm4°^ = -1.2394n/eF, (3.69) 

d^^hS^^^ = -0.6290n/eF, (3.70) 

dhh£^^^ = -1.2394n/ej., (3.71) 

ei^/^^ = -0.0509nei., (3.72) 

e(i/^) = -0.477nei., (3.73) 

a^e(^i/JV) ^ _o.513n, (3.74) 

a;,4VJV) ^ o.452n, (3.75) 

^Aei^/^^ = -0.379n, (3.76) 

a^e(i/JV) = _o.948n. (3.77) 



field vanishes, (2) the normal state by definition does not have 
any superfluid correlations, i.e., its ground-state energy is in- 
dependent of A, (3) the ground-state energy of the superfluid 
state is a minimum with respect to A. We then obtain 



(3.60) 



(3.61) 
(3.62) 
(3.63) 



At this point a few remarks are in order. Some of these 
matrix elements can be related to each other Indeed, using 
scaling arguments akin to those presented in Sec. lIII A~3l one 

can show, for instance, ei^^^' = l/icJ^ei^^^'' + IAc^a^s^^^'', 

as well as e"n^^^ — \[i'dfj,e"n^^ + ^hdh£n^^\ Furthermore, 
we observe that the corrections Sfi and 5h do not involve any 
matrix elements containing derivatives with respect to A. This 
merely reflects the fact that A is a variational parameter and 
that one could have reduced the solution for Sfi and Sh to a 
set of two coupled equations using Eqs. ( |3.50t . In addition, 
we point out that the ratio A//i at /i = hc2 is equal to A//i at 
h — 0, a result directly related to the vanishing susceptibility 
of a gapped superfluid state. 

The polarization at the upper-critical field can be written in 
a 1 /A^-expansion as 



Pc.^^= + 4<5P, (3.78) 



where is the result of the N ^ oo calculation and 



ff, JO) ^ (0) ^ (0) 

dAt^ef^ C^AAE 




Solving the matrix equation (13.601 1 yields 



5h 
(5A 



I 
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6P = Pi°) 



dh£ 



(l/N) 



(0) 



(l/N) 



dnhS 



(0) 



(0)' 



5 .(0) 



(0) 



{or 



dhSn 



(0) 



5 .(0) 



5h 



(3.79) 



At unitarity, we find 

0.6930 + 0.087/7V + 0{1/N^), 
0.9326 - 0.631/Ar + 0{1/N^), 
0.8586 - 0.458/7V + 0{1/N^), 
0.9978 - 0.242/7V + 0{1/N^). 



Substituting previous equations into Eq. ( 13.311 1 yields, 



hc2/eF 



(3.80) 
(3.81) 
(3.82) 
(3.83) 



C. Away from unitarity 



Although so far we have focused on the problem at the uni- 
tarity limit, the same large-A^ expansion formalism can be 
used to calculate the properties of a resonant Fermi gas away 
from the unitarity point. In this section we compute such 
correction at zero temperature. 

The result of Eq. (13.30b holds for any value of the scattering 
length as and therefore can be used away from unitarity. The 
matrix elements of the left hand side of Eq. ( |3.30t can be cal- 
culated by taking appropriate derivatives of Eq. ( 13. 2t . After 
straightforward algebraic manipulation, one obtains 

















= l5[Xo) 



(2m) 



3/2 



27r2 
(2m)3/2 
27r2 



Ai/2 
Ai/2 



)3/2 



27r2 



(3.84) 

(3.85) 
(3.86) 

(3.87) 



where, following Marini, et al. — , we introduced dimension- 
less function l5{xo) and Ie{xo) 



h{xo) ^ dx 



x% 
El 

2 



X 



An 2m 



(3.88) 
(3.89) 



and dimensionless parameters 

L^^^x^-xo , E^ = ^^^e. + 1. (3.90) 

The variables /io and Ao are the solutions to the gap and 
particle equations in the N ^ oo limit. For numerical pur- 
poses, it is useful to note that Eqs. (13.881 1 and (13.89b can be 
expressed in terms of linear combinations of complete elliptic 
integrals^. 



2e 



3/2 



1 



^^'~3nAl/'m^o) + Iiixo) 



for the leading order (in l/N) corrections 5 A and (5/i to the 
saddle-point (mean-field) result given by Eqs. ( 13.31 ) and (13.41) 

at T ^ 0. 

Solving Eq. ( 13. 91b numerically yields results displayed 
FiglHas a function of {kpag)^^ . 



IV. DISCUSSION AND COMPARISON TO OTHER WORK 

In this work we have studied a two-component resonantly- 
interacting Fermi gas in the unitary regime by generalizing 
the system to a 2A^-component gas with N flavors of spin- 1/2 
fermions, computing a variety of thermodynamic quantities 
in a systematic expansion in a small parameter 1 /N, with our 
main results presented in Sec. IIBI Although our results are 
strictly only quantitatively valid to leading-order in 1/iV, to 
compare to experiments and other calculations, in this section 
we boldly take the iV ^ 1 limit. 

The study of the unitarity point can be seen as a benchmark 
for many-body theories. In this section, we collect the various 
estimates to thermodynamic quantities at the unitarity point 
from the literature and compare to our results. Starting with 
BCS expressions, a mean-field estimate of ^ = 0.59 can be 
obtained. Since this approach is based on a variational wave- 
function, it provides a strict upper bound. Pade approximant 
techniques have also been applied to a Fermi gas expansion in 
terms of kpas to extract an estimate for ^ = 0.326'^'^^. Perali, 
et alP have introduced a diagrammatic method based on the 
t-matrix extended to the superfluid phase finding, at unitarity, 
^ = 0.455. Most recently, Haussmann, et al?'^ proposed a 
self-consistent and conserving theory to study the crossover 
and found ^ = 0.36. An epsilon expansion of the unitarity 
gas, based on an expansion from the dimension d ~ A — e 
has been performed to second loop order™. Depending on 
the Borel-Pade extrapolation schemes used, the results range 
from ^ = 0.30 to ^ = 0.37. 

At present time, the best estimates for ^ coming from fixed 
node Green's function Monte Carlo calculations''''^" yield 
i w 0.44 ± 0.01 (by Carlson, et al}'^) and i w 0.42 ± 0.01 
by Astrakharchik, et al^. The nodal method is based on a 
variational approach, thus also giving a strict upper-bound for 

e 

There has also been a large experimental effort aimed at 
extracting the parameter ^. Experimental measurements of 
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TABLE I: Recent experimental results for f compared with calcu- 
lated values (compilation taken from Ref. |23 with a few additional 
results) 











Gehm et ai— 


0.74(7) 




Bartenstein et al— 




Experimental 


Bourdel (2004) et al^ 


0.36(15) 


results 


Duke (2005)^' 


0.51(4) 




Partridge et al^ 


0.46(5) 




Regal et al— 


0.38(7) 




Stewart et al.^- 


'-'•^"-0.12 




Mean-Field 


i < 0.5906 




Astrakharchik et al. 


0.42(1) 




Carlson et ai— 


0.44(1) 


Calculated 


Perali et al.— 


0.455 


values 


Fade approximation^^ 


0.33 




Haussmann et ali^ 


0.36 




Epsilon expansiori^ 


0.30 to 0.37 




This work 


0.28 



TABLE II: Universal ratio Acxc/ef at unitarity: Comparison be- 
tween numerical and theoretical approaches. 





Acxc/eF 


Mean-Field 


0.6864 


Carlson et al. (Monte Carlo)^ 


0.54 


Haussmann et al.'^ 


0.46 


Nishida et a/.— (epsilon expansion) 


0.60 


This work 


0.49 



the expansion of ^Li in the unitarity limit from a harmonic 
trap have determined ^ (often quoted as /? = ^ — 1). Exper- 
imental results for ^ were first obtained by Ghem, et al^, 
i = 0.74(7) , Bai-tenstein, et al.'^ S, = 0.32l}i] , and 
Bourdel, et al}, ^ = 0.36(15) . Most recently, the Duke 
group^' obtained ^ = 0.51(4) , whereas the Rice Group'*' 
found ^ = 0.46(5) . Experiments on ^"K by the Boulder 
group ^'^ yielded ^ = 0.46j:°:'j'2- Our prediction for ^ = 0.28 
at = 1 is qualitatively consistent with these results, re- 
flecting the smallness of fluctuations around the saddle-point 
mean-field solution at T = 0. 

Another important quantity is the single particle excitation 
gap at unitarity. In the Green's function Monte Carlo calcula- 
tions of Refs. 19 53, an estimate for the spectroscopic energy 
gap Acxc was determined from the odd-even staggering of the 
total energy as the number of particles is increased and found 
Acxc — 0.54ei?. It is not a priori clear whether this gap is 
generally identical to the expectation value A = |A(7/'|-0i) |; 
in the Appendix we show that they are equal to zeroth order 
but differ at the 0{\/N) computed here. Our prediction for 
Aoxc/fiF, evaluated at iV = 1, is in good agreement with the 
Monte Carlo calculations as seen in table HIl 

The critical temperature at unitarity has also been inves- 
tigated. The 1/A^ calculations shows relatively large cor- 
rection to the ^ oo result, so that the ^ 1 limit 
of Eq. ( 11.81 ) is clearly not sensible quantitatively but can be 



qualitatively interpreted as predicting a large negative cor- 
rection to the mean-field result for Tc. This strong correc- 
tion to the saddle-point approximation is expected and re- 
flects the fact that the mean-field solution neglects the effect 
of bound pairs and describes the normal state as consisting 
of free fermions^. Quantum Monte calculations performed 
by Bulgac, et ali^ have indeed determined a relatively low 
critical temperature (Tc = 0.23(2)6^), and Burovski, et al.^ 
have arrived to Tc = 0.152(7)ei? using diagrammatic deter- 
minant Monte Carlo. We note that this result is below the 
BEC limit, Tbec ~ 0.218eF, in contrast to the earlier work 
of Nozieres and Schmitt-Rink^'. Burovski, et al. also give 
/i = 0.493(14)eF and e = 0.31(l)eF at T^. Finally, we point 
out a recent 1/A^ calculation by Nikolic and Sachdev^^, that 
gives a result for 1 /Tc, that is consistent with our prediction 
ioxTc. 

Now we turn to the polarized Fermi gas. We have deter- 
mined the 1/iV correction to the upper critical field Pc2- Be- 
sides the mean-field value, we do not know of any other the- 
oretical estimates for Pc2. This is partly due to the difficulty 
of Monte Carlo calculations to tackle the sign problem for a 
polarized Fermi gas. However, a recent experiment on ^Li has 
found Pc2 to be below its mean-field value estimate, and given 
by Pc2 = 0.70 This is consistent with the calculation of the 
1/7V correction that shows that Pc2 decreases with decreasing 
N going from its mean field value P^^°° = 0.933. The naive 
substitution iV = 1 in Eq. [IITj vields Pj^=i = 0.302, which 
underestimates the experimental result. 



V. SUMMARY 

To summarize, we have studied a resonant Fermi gas inter- 
acting via short-range attractive interactions by generalizing 
the model to N fermion flavors and employing a theoretical 
method that is perturbative in 1/N . The 1/iV expansion pro- 
vides a systematic scheme for quantitatively determining cor- 
rections to the standard BEC-BCS mean-field theory of inter- 
acting Fermi gases. 

Although our primary goal was the computation of various 
quantities to leading order in 1/A^ in the vicinity of the uni- 
tarity point (where universality holds), we also computed the 
zero-temperature gap and chemical potential away from the 
unitary point. Clearly, future work must be done to gener- 
alize our other results for arbitrary coupling (i.e., Feshbach 
resonance detuning) as well as to compute the next 0(l/iV^) 
term. 

Note Added: We thank Y. Nishida and D. T. Son for point- 
ing out an error in an earlier version of this manuscript. We 
learned of an independent, and closely related 1 /A^-expansion 
study by R Nikolic and S. Sachdev^^. Where there is overlap, 
our results are in agreement with theirs. 

We gratefully acknowledge useful discussions with V. Gu- 
rarie and A. Lamacraft, as well as support from the NSF Grant 
DMR-0321848. 



12 



APPENDIX A: PROPERTIES OF GREEN'S FUNCTIONS 

In this Appendix we determine properties of the single- 
particle excitation spectrum. In particular, we calculate the 
single-particle excitation gap Aoxc to subleading order in 
1/A^). Although Acxc (which is more commonly measured in 
experiments) is equal to the order parameter A within mean 
field theory (i.e., at iV — > oo), beyond mean-field theory they 
are slightly different. 

The gap in the single-particle excitation spectrum is deter- 
mined as the lowest energy state of the single-particle Green's 
function, defined as 



where the Green's function matrix is given by 

b{x) 



-9r + h+ /i 



— Or — 77 U J 



(A.8) 



Using the decomposition given in (Eq. 12.131) . we find to sub- 
leading order 



{{G{k)) 



sib] J 



where 



(A.9) 



(A.l) 



(A. 10) 



where without loss of generality, we selected the flavor i = 1. 
To calculate this expectation value within our formalism, it is 
useful to introduce source fields in the action 



In particular, we find for the self energy 

1 



a 



source 



^ r(i r _ _ 

J2 dT d^r [$,(x)*,(x) + $,(a;)4',(a;)] , 
,;=l Jo J 

(A.2) 

where we introduced a Nambu representation for fermions, 

i.e. ^i{x) = (V'i-f (x), and where ^i{x) = 
(77j|(x),77ij^(x)) are the source fields. From the auxiliary 
fields, we derive 



N ^ \A{k - k')\^ - \B{k - fc')p oj'^ + El, 

A{k' - k){iJ - i]^,) -B{k~k')A* 
-B*{k- k')A -A{k - k'){iLu' + ^k' )^ 

(A. 11) 



By solving the equation det^G^k, uj))~'^]iuj^uj — 0, we 
obtain the dispersion relation. At leading order, i.e. N 00, 
this is given by 



92 lnZ[?7, rj] 



drii,a'iO)drj,Jx) 



(A.3) 



7^ — 77—0 



Using the formula 

.(2,)£)^^(2,)g-2:.J*i(2:)Jl^(^)*.(^)+*i(2:)*.(2:)+*>(a:)*>(x)] 

(A.4) 



we find the partition function is given by 



Since the chemical potential is positive, it is possible to find a 
wavevector such that the condition = is fulfilled, giving, 
at leading order, that the lowest energy excitation is equal to 

the order parameter i.e. Acxc = A^"'. At subleading order, 
the additional contribution due to the self energy breaks this 
equality. In this case, we determine the excitation gap by ex- 
panding the dispersion relation around its minimum, namely 
Ek = M- Straightforward manipulations give 



where 



] = Z-^ J D6(x)m*(a;)e^^[''l-^=[*1, (A.5) 



Acxc = ^ + l (^11 + ^22 - 2S12) 



(A. 13) 



^ dr I d^r<^>,{x)G{x)^,{x), (A.6) 



where A is the order parameter calculated at order 1 /N. Nu- 
merical evaluation of the self energy gives E11+E22 — 2S12 = 
— 0.067eF, yielding 



and S[b] is given by Eq. ( 12.101 ). Using, Eq. ( IA.3I ), we arrive to 

g{x)^{[G{x)],,)s[b], (A.7) 



Acxc/ej^ = 0.6864 - 0.196/7V + 0(1/^^2), 



the final result. 



(A. 14) 
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